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The Nambu-Jona-Lasinio model reduced to 2+1 dimensions has two different path integral for- 
mulations: at finite chemical potential one formulation has a severe sign problem similar to that 
found in QCD, while the other does not. At large N, where N is the number of flavors, one can 
compute the probability distributions of fermion correlators analytically in both formulations. In 
the former case one finds a broad distribution with small mean; in the latter one finds a heavy tailed 
positive distribution amenable to the cumulant expansion techniques developed in earlier work. We 
speculate on the implications of this model for QCD. 
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QCD has been around for forty years and yet its application to the properties of bulk matter at nuclear density 
has proven to be stubbornly intractable. An inherently nonperturbative problem, the only reliable tool available for 
solving QCD in this energy regime is Monte Carlo evaluation of lattice QCD. Unfortunately, lattice methods suffer 
from a severe sign problem in the grand canonical formulation which renders Monte Carlo techniques useless. In a 
canonical formulation, despite recent impressive progress in studying light nuclei and hypernuclei p], there remain 
severe problems with signal-to-noise; these problems are clearly related to the sign problem. We continue here the 
research program outlined in Ref. [2H1] where we study the probability distributions of correlators and argue that the 
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I. INTRODUCTION 
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origins of the sign or noise problem lie in the dynamics and spectrum of the theory, and that it is not especially useful 
to think of it as a "fermion" sign problem. In the particular case of QCD, the severity of the sign problem is closely 
related to the phenomenon of chiral symmetry breaking and the consequential existence of a light pion. 

In this paper we elucidate the connection between the QCD sign problem and chiral symmetry breaking by studying 
a simpler theory - the Nambu-Jona-Lasinio (NJL) model [5] in three dimensions. This formulation of the NJL 
model with N flavors is of particular interest because it is soluble in large N, because it exhibits chiral symmetry 
\Q breaking without the complication of confinement, and because it has two complementary but equivalent path integral 
formulations - one of which looks very QCD-like and has a severe sign problem, while the other resembles more closely 
the chiral quark model [6] and has no sign problem. 

We begin by reviewing some features of the sign problem in QCD, then turn to our analysis of the NJL model. In 
particular we are able to analytically compute probability distributions for correlator measurements and show how the 
noise spectrum in these measurements is related to the presence or absence of a sign problem in the grand canonical 
formulation. Our results should be directly applicable to suitable lattice formulations of this theory, but all of our 
^ analysis is analytic and in continuum Euclidian spacetime |22) . 

At the end we return to QCD and speculate on the implications of our analysis for finding a path integral formulation 
of QCD which would allow for numerical study of bulk matter. 



II. THE SIGN/NOISE PROBLEMS IN QCD AND THE UNIQUE ROLE OF THE PION 

Numerical simulations of lattice QCD involve an approximation of the Feynman path integral in Euclidian spacetime. 
This requires a Monte Carlo evaluation of the averages of operators of interest - such as hadron correlation functions 
- over an ensemble of background gauge field variables generated with probability distribution 

V{U) oc e- s ^A[U] (1) 
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where U is a link variable for the gluon degrees of freedom, S is a suitable discretization of the Yang-Mills action, 
and A is the fermion determinant; both S and A are functionals of the link variable U. This program has been very 
successful for determining properties of the QCD spectrum in the vacuum, but runs into an obstacle when trying to 
simulate matter at finite density in a grand canonical ensemble. The fermion determinant A is a discretized version 
of det(0 — to + fJ/ji), where is the covariant derivative, m is the quark mass, n is the chemical potential for 
quark number, and the gamma matrices 71, . . . , 74 are all hermitian. At nonzero /x the fermion operator is in general 
complex, since fxji is hermitian while Ip is anti-hermitian, and the two do not commute. As a result the expression 
in eq. is not a suitable probability measure for the gauge field configurations; this problem has nothing to do 
with discretization of the lattice action, since it is a property of the continuum functional. Therefore, for the rest of 
this paper we will only discuss the sign problem in the continuum, as the discretization of spacetime is not directly 
relevant (although a poor choice of discretization can in principle introduce additional sign problems not present in 
the continuum theory). 

At nonzero chemical potential one can write A(A fl ,pi) = |A(-A M , p)\e l6 , and one finds that the phase starts to 
fluctuate wildly with changes of the gauge field for fi > m n /2 at low temperature [7]. This behavior seems precocious, 
since at T = 0, nonzero p can have no effect on the free energy until p > totv/3 > m„/2, where tojv is the nucleon 
mass. The phenomenon was clearly explained by Splittorff and Verbaarschot [HI [S] (SV) who noted that for two 
degenerate flavors, ^(A^, //)| correctly describes the fermion determinant with a chemical potential p for isospin — 
namely +p for the up quark and — fi for the down quark. Such a system will exhibit Bose-Einstein condensation of 
pions, with free energy becoming rapidly negative for p > m n /2, as the pion is the lightest state carrying isospin. From 
chiral perturbation theory they derive the estimate (in a continuum Euclidian spacetime volume V, for p > m^/2) 

2fxFn/ml e _2VF> 2 (l-rr4/V) 2 ( 2 ) 
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where F w is the pion decay constant. This result shows that the phase 9[A) is responsible for cancellations in the 
integral that lead to this ratio going to zero exponentially with the volume for p > m v /2. The sign problem can be 
thought of as being necessary to keep pions out of the ensemble for nonzero baryon number. 

A complementary obstacle is seen when studying the baryon spectrum in a canonical formulation of lattice QCD, 
at zero chemical potential. A typical approach to compute the mass Mb of the ground state with baryon number B is 
to use Monte Carlo methods to measure the correlation function Cb(t) — (O(0, r; A^)), where O = G(0, r; A M ) 3S , G 
being the quark propagator from time t = to t = r, with color, flavor, and spin indices appropriately contracted to 
produce a state with the desired quantum numbers. This correlations function must behave as Cb(t) oc exp(— Mbt) 
at large r and one can therefore extract Mb as the limit of Mb = — In Cb(t)/t as r — > 00. This procedure seems 
straightforward and innocent of any sign problem since we have set p = and the fermion determinant is real and 
positive. However, one finds that a Monte Carlo measurement of Cg(r) is very noisy at late time. A particularly 
simple explanation was given by Lepage [10) : he pointed out that the variance in the measurement of Cb is given by 
the average a 2 = (G(0, r; A tl ) 3B G 1 > (0, r; A p ) 3B ), corresponding to the propagation of 3B quarks and 2>B anti-quarks, 
with a ground state consisting of 3i? pions. Therefore one would expect a 2 cx exp(— SBiti^t) at late time, with the 
signal to noise ratio for the Monte Carlo measurement of Cb being proportional to exp(— 3t(Mb/3 — Brn-^/2)), which 
vanishes exponentially fast since the mass of B pions is less than 2/3 the mass of the i3-nucleon ground state. Thus 
eventually the signal will always be overwhelmed by noise, and again it is due to the pion being lighter per constituent 
quark than the nucleon, just as we saw in the grand canonical example. If one works at fixed baryon density p with 
B = Vp, where V is the volume, one sees that the noise problem grows exponentially with the volume, as one would 
expect from the problem encountered in the grand canonical approach, eq. p]). 

Savage has extended the analysis to look at higher moments of the distribution function for Cb [23] . Even moments 
all involve equal numbers of quarks and anti-quarks, and therefore fall off exponentially with a rate determined by 
the pion mass. However odd moments involve expectations of operators with net baryon number, and fall off more 
quickly, relative to a — exponentially with the nucleon mass. With odd moments vanishing faster with r than 
even moments, one can conclude that the probability distribution for Cb evolves exponentially fast at late time to a 
symmetric distribution with vanishing mean, a manifestation of the sign problem which in either approach is due to 
the existence of the light pion. 

There have been a number of microscopic analyses of the QCD sign problem (for example, Ref. [5]), as well as 
general proposals to modify the conventional approach to simulating the Feynman path integrals (such as the meron 
cluster approach, Ref. Here we will instead pursue further the connection between the sign problem and chiral 

symmetry breaking, the origin of the pion's small mass. 
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III. THE LARGE-TV NJL MODEL IN THREE DIMENSIONS 



We consider a modified version of the original NJL model 5 , with TV flavors and reduced from four to three 
Euclidian dimensions: 

C = N(^ a (0-m)lP a - | [(Ma) 2 + &aH5lpa) 2 ]) ■ (3) 

Our convention is that a,b . . . are flavor indices summed over 1, . . . , TV, the 3D coordinate indices are ., summed 
over 1,2,3, while Greek indices fi, v... are summed over 4D coordinates 1, ... ,4. Thus $ = jidi, for example. The 
gamma matrices are the usual 4x4 matrices used in 4D, not the 2x2 matrices appropriate for 3D, and so the 
Lagrangian represents 27V flavors of 3D Dirac fermions. However, we will generally use 4D nomenclature — in 
particular, in the limit m this theory has a U(l) chiral symmetry in 4D, which becomes a flavor symmetry in 
3D; despite the fact that there is no notion of chiral symmetry in 3D, we will continue to refer to this global U(l) 
symmetry as a chiral symmetry. 

In 3D this model still exhibits the key feature of the original NJL model, the phenomenon of spontaneous chiral 
symmetry breaking which gives rise to a Goldstone boson we will call the pion. We will analyze this theory in a 1 /TV 
expansion using techniques similar to those used in the Gross-Neveu model, which is formulated in 2D [23]. This 
theory has been reviewed in jTSJ [13] and has been studied numerically in [14] . 



A. The ct/tt formulation 



The conventional treatment of the theory eq. ([3]) is to introduce auxiliary fields a and 7r to obtain a bilinear fermion 
action, 

C = N Q- (a 2 + 7T 2 ) + V> [$ - m + a + in l5 ] ^ . (4) 

In this formulation, the a and it fields are singlets under the SU (TV) flavor symmetry, while <fi = (cr+iir) / v2 transforms 
linearly under the chiral U(l) symmetry. With the above normalization, TV-counting is simple: every vertex and every 
fermion loop brings a factor of TV, every propagator a factor of 1/TV. Loops which include scalar propagators do not 
give a factor of TV, since the mesons do not carry the TV flavors. 



1. No sign problem 



An interesting feature of this theory is that at finite chemical potential the fermion determinant is real, and for 
even TV, it is positive. To see this, note that we can define a real symmetric charge conjugation matrix C that satisfies 
C 2 = 1, Cj^C — 7* for i = 1,2,3 and C75C = 75 . For example, we can take 7* = a\ x at, 75 = 03 x 1 with 
C = cr% x 1T2 • Then the fermion operator for a single flavor in the grand canonical formulation satisfies D* — CDC, 
where D = ($ — m + a + wr7s + /U71). It follows that complex eigenvalues of D come in conjugate pairs, while real 
eigenvalues can be unpaired. Thus det D for the TV-flavor theory eq. Q equals a real number raised to the power 
TV and is positive for even TV. This implies that there is no sign problem obstacle to simulating this theory at finite 
density using Monte Carlo methods [T5] . 



2. Chiral symmetry breaking 



Integrating out the fermions gives rise to the effective action S e s = NS(<j,tt), where 

" 1 



S[a,*] = J d 3 



2/T 



(cr(x) 2 + 7r(x) 2 ) — Tr ln($ — m + &(%) + i^(x)j^) 



(5) 



The large-TV expansion is equivalent to the semiclassical expansion, and the vacuum is characterized by the clas- 
sical solution that minimizes the action S. We can readily compute S(a, 0) with a — constant using dimensional 
regularization and minimal subtraction (MS) |25j . Up to an irrelevant additive constant we find 



S(ct,0) = VT 
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(— m + er) 3 + m 3 
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(6) 
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where we have put the system in a box of spacetime volume VT. Defining the chiral symmetry breaking minimum to 
be at (a) = f when to = 0, we obtain 



dS(a,0) 



da 
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9 = - 



m=0 



f 



(7) 



With this value for g (which is renormalization scheme dependent), we find for nonzero quark mass to > 

" <r 2 (2<7 - 3/) + 6ma(m - a) 



S(a,0) = VT 



6ir 
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with minimum at 



1 + v/f+4?) + 2r\ 



m 



f ' 
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which shows how the chiral symmetry breaking minimum depends on the explicit quark mass. 

Next we expand the effective action S to second order in spacetime dependent fluctuations of the a and tt fields, 



S(a,n) = S((a),0) + X -V 2 T 2 



d 3 k 
(2tt) 3 



[D a (k)5u(k)5a(-k) + D n (k)5n(k)6ir(-k)} + ... 



(10) 



where D a and D v are readily computed in the MS scheme from one-loop diagrams in the background (a). The 
constituent quark mass, M, is given by 



f r 

M = {a) - to = - 1 + y/l + 4?7 



(11) 



Since there is no confinement in this theory, M is the mass of the lightest fermionic excitation and to is the current 
quark mass. With the fcrmion propagator, G(p) = + M) -1 , we find, 



ArM = - g + 
1 



D„(k) = 



fJ 



d 3 q 
(2tt)3 

d 3 q 

(2tt) 3 



Tr [G(q + k/2)G(q-k/2)} = - 
Tr [ l5 G(q + k/2) l5 G(q - k/2)\ 



f (AM 2 + k 2 ) cot" 1 (2M/k) + 2Mk 
7r 2irk 
_ 2M - 2/ + fccot- 1 (2M/k) 
2tT~ ' 



(12) 



To leading order in 1/N, these functions determine the a and tt dispersion relations, and their masses are defined by 
the location of zeros in D a and D^, respectively. One finds m 2 = (2/) 2 in the chiral limit. However, the pole sits at 
the beginning of the two-fermion branch cut, and for to > we find m a > 2M and so the a field is unstable. Only at 
subleading order in 1/N would one see the branch cut appear in D a for a — > irir decay, as that entails an additional 
quark loop. A chiral expansion for the pion mass yields 
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and near the pion pole, k 2 = —m 2 , one finds for the pion propagator 
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ND^k) Nk 2 + m 2 
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(13) 
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The positivity of m 2 and m 2 for to > shows that we found the correct (e.g. stable) vacuum. 



B. The A/V formulation 



An alternative formulation of the theory follows if one first uses the Fierz identity 

[SijSki - (75)ij(75)w] = \ [{iphilnhj ~ {7^5)il(^l5) kj] 



(15) 
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to rearrange the four-fermion interaction in eq. (|3j) before introducing auxiliary fields. The Lagrangian becomes 

<i 

where $ = jidi with i summed over i = 1,2, 3, while the 7^ matrices are summed over \i = 1, . . . , 4. This formulation 
invites introduction of N x N matrix valued vector and axial vector auxiliary fields V and A, giving the equivalent 
theory 

C = N Q Tr (V^ + A^) +^[0-m + if + A ■ (17) 

iV-counting in this theory is different from the a /ir formulation, since the A and V meson fields are N x N matrices. 
In fact, the N counting here is identical to that of large- N QCD, and it is convenient to employ 't Hooft's double-line 
notation for the mesons. As in QCD, the order of a graph without external legs is given by N x , where x 1S the Euler 
characteristic of the surface defined by the graph, and so to leading order one only need consider planar graphs with 
a minimal number of closed fermion loops. However this class of graphs is much simpler than in QCD, since the A 
and V mesons have no cubic or quartic interactions, unlike gluons. 



1. A QCD-like sign problem 



In this A/V formulation, the fermion matrix at finite chemical potential is given by D(/j.) = 
\0 — m + + -^75 + M7i] i which is similar in structure to the QCD Dirac matrix with nonzero chemical poten- 
tial [i (with the N flavors playing the role of color) and its determinant is similarly complex. In fact, as in QCD we 
can make an SV argument [HI [S] about the phase of the fermion determinant by considering two degenerate families 
(e.g. 2N fermions) so that the chiral symmetry is enlarged from Ul(1) x Ur(1) to Ul(2) x Ur(2). The fermion 
determinant in the case of a quark number chemical potential is (det _D(/i)) 2 while with an isospin chemical potential 
it is I det D(p)\ 2 , the difference between the two being the phase e 2lB . In the latter case there is a transition to a 
pion-condensed state at fi = m^/2 just as in QCD, and so the SV argument leads a similar formula as in eq. pi). 

It is remarkable that a theory with a sign problem so similar to that of QCD is known to have a formulation which 
has no sign problem, the a fix formulation of the previous section. 



2. Chiral symmetry breaking 



One cannot see chiral symmetry breaking in the A/V formulation in a mean field formalism, as there is no fun- 
damental field with the right quantum numbers to play the role of the ijjtp condensate — again, as in QCD. Instead 
one has to consider the Schwinger-Dyson equation for the fermion propagator, which is exact in leading order in 1/N 
since vertex and meson propagator corrections occur at higher order. We take the full canonically normalized fermion 
propagator to equal 



which satisfies the integral equation 
-ifZ{p)+M{p) = - 

— —ifl — m + 4g 
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d 3 k 
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(18) 
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-ifkZ(k) + M{k) 



7^75 



(2ir) 3 k 2 Z(k) 2 + M(k) 2 



(19) 



From this one finds (using the MS renormalization scheme as before) that Z(k) — 1 and M(k) = M is a constant 
satisfying 



M 



-m + 4g 
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(2tt) 3 k 2 + M 2 



gM 2 



which has the solution 
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FIG. 1: A graphical representation of the integral equation eg. (23) for the four-point correlator in the A/V formulation for an 
incoming fermion/anti-fermion pair of one flavor and an outgoing pair of another. Dirac indices are labeled. 



Using the renormalization condition that the dynamical fermion mass equals / in the the chiral limit, m = 0, gives 
g = —71"//, as in eq. ^ and the solution for M 



M = 



f 



(l + V 1 + V) V 



m 

1 ' 



(22) 



which is the same as the value derived in the a fix formulation for the dynamical fermion mass, eq. (11 1 



In the A/V formulation the a and n mesons appear as a fermion and anti-fermion pair bound together by strong 
and V„ vector meson exchange, much the same way as mesons arise in large-iV QCD. To see them we compute 
the connected four-point function Ai for a fermion/anti-fermion pair of flavors £,k, each with 3-momcntum k/2, to 
scatter into a fermion/anti-fermion pair of flavors i,j with 3-momenta k/2 ± p. To leading order in l/N, Ai obeys 
the integral equation, shown graphically in Fig. [T] 
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2 J (27T, 
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(23) 



where G{p) is the free fermion propagator with dynamical mass M , and we used a Fierz identity to replace the vector 
and axial vector gamma matrices by scalar and pseudoscalar. The solution to this equation is 



Mij-ki{k) = -N 



Sjjhi {ilb)ij{ils)ki 
D„{k) D„{k) 



= -N 2 (SijSkiG^k) + (*7B)«(»78)fciG w (ifc)) 



(24) 



where we dropped the label p (as our solution is independent of p), D a lx are given in eq. (12), and G a ^ = N/D a ^ are 



the full meson propagators. Thus we see that (up to a sign) the interaction between valence fermion and anti-fermion 
via f-channel exchange of A^ and mesons is exactly equivalent to an annihilation diagram in the a/n formulation, 
with a single meson in the s-channel; similarly, the interaction between two valence fermions or two valence anti- 
fermions in the A/V theory is equivalent to a single meson exchange in the it-channel in the a/n theory (Fig. [2]). This 
equivalence will allow us to use the simpler a/n theory to calculate the cumulants for the A/V theory. 



IV. THE PROBABILITY DISTRIBUTION OF THE FERMION PROPAGATOR 



If X[4>] is a functional of a stochastic field 4> corresponding to an observable — such as a correlation function — we 
define the normalized probability density function (pdf ) for X to be the path integral 

V(x) = M [ [d0] e~ s ^ s(X[<l>] - a:) . 



(25) 

where we assume S[<p] is real. If one were to sample an ensemble of <j> configurations according to the distribution 
e -S[4>]^ values of X[<fi] would be distributed according to V, making its relevance to Monte Carlo simulations 
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FIG. 2: In the A/V formulation at leading order in 1/N, interaction between a valence fermion/anti-fermion pair (M) or a 
valence fermion pair (M' ) is equivalent to s- or u-channel exchange respectively of single a and it mesons in the a/n formulation, 
where the fermions have mass M arising from nonzero (a) . 



evident. For an accurate estimate of (X) from a reasonable number of samples, one would like V(x) to be sharply 
peaked around its mean. However, one might find a very broad distribution centered about a mean close to zero, 
making an accurate estimate of the mean, without a huge number of samples, very difficult, such as is the case with 
baryon propagators in QCD. Alternatively, one might find a heavy-tailed distribution for which very rare events make 
a significant contribution to the mean, resulting in very noisy and often misleading estimates of (X) from a finite 
sample. This latter situation is indicative of what is called "an overlap problem", which occurs when e~ s ^ is peaked 
far from the field configurations which provide support for nonzero X((f>). With some knowledge about the nature of 
the tail of the distribution, it may be possible to use statistical methods to greatly improve the determination of (X) 
0- 



The pdf given in eq. ( 25 ) is a difficult quantity to analyze using field theoretic methods because of the singular 



nature of the delta function; instead, we consider the characteristic function (cf) $(s), which is just the Fourier 
transform of the pdf: 



(26) 



This is a useful formulation because W(s) = — ln$(s) is on the one hand given by the connected Feynman diagrams 
[26| of the modified action S[<p] —isX[</)], while, up to an irrelevant additive constant, it is also the generating function 
for the cumulants of X: 



W(s) 



(i S y 



n=l 



(27) 



Here, K n is the n th cumulant, with n\ = (X) = fi being the mean of V{x), K2 = ((X 2 ) — (X) 2 ) = a 2 being the 
variance, etc. 

This procedure has to be modified slightly when dealing with a complex observable, where we replace eq. ( 26 1 by 

s) = e - w ^ = TV J [&<t>] e -sm+Ksxm + -sxm , (28 ) 

where s is now complex and the bar indicates complex conjugation. Now, In W(s, s) has a double expansion in both 
s and s, 



W(s,s) = - J2 



(is) m {is) r 
ml n! 



with K n _ rn — K m ,n- For example, 
«i,0 = ( x ) , 



n hl = (\X\ 2 )-\(X)\ 2 , K2 , = (X 2 )-(X)' 



(29) 



(30) 



In what follows we show how to compute the cumulants of the correlation function of a single fermion with zero 
2-momentum for a time extent r. This is the sort of measurement one would perform in a Monte Carlo simulation in 
order to determine the mass of the fermion. In particular we take for X[<f>] 



Y r = hi 



0,t = r/2|Trr 



i 



*r = ^(p = 0,< = T/2|Trr 



— m + a 
1 
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+ if + 4 75 
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a /it formulation, (31a) 
A/V formulation, (31b) 



where T is some Dirac matrix of our choosing. Our choice to look at the log of the propagator in the a/n formulation 
makes the later analysis more simple. Note that in the definitions of our observable Xr and Yr we use canonically 
normalized fermion propagators, without the factor of 1/N. Measuring the expectation value of this correlator is a 
procedure for determining the mass to/ of the lightest fermion state through the formula 

lim --ln(X r )=lim -- ln(p7 r > = m f , (32) 
provided that V does not project out this state. Of course, we already know analytically that to/ = M, with M 



given in eq. (11), but by calculating the cumulants for this observable we will establish how difficult it would be 
to determine the fermion mass by numerical Monte Carlo methods using the two formulations of the NJL model, 
and why. In particular we will show that in the QCD-like A/V formulation, the pdf for Xr, at late time, looks as 
one would expect from the Lepage-Savage picture: a broad distribution that is nearly symmetric about zero with 
an exponentially small mean. In contrast, the pdf for the physically equivalent a/n formulation looks heavy-tailed 
and close to lognormal at late time. Thus, a Monte Carlo study of this theory, without a sign problem, still faces an 
overlap problem and significant numerical challenges, but is perhaps amenable to a cumulant expansion as introduced 
for a similar system in Ref. [2H1] ■ This theory has also been successfully investigated recently using the "fermion bag 
approach" [16HI8] . 



A. Noise distribution in the A/V formulation and the Lepage-Savage analysis 

1. The tiifi and cumulants 

We begin by computing the cumulants for measurements of the fermion propagator Xr in the A/V formulation; 
since this observable is complex, we use the formalism in eqs. ( 28||30 1. From the above discussion, the K mn = K nm 
cumulant for Xr is given by the sum of connected graphs with m copies of Xr and n copies of its complex conjugate, 
Xr. We will refer to these as valence fermion and valence anti-fermion propagators, respectively; at leading order in 
1/N there are no sea quark loop contributions. With our definition of X, there are no factors of 1/N from valence 
fermion propagators, nor factors of N from meson coupling to valence fermions, nor is annihilation of a valence fermion 
with a valence anti-fermion allowed. The computation of n m ,n involves graphs with net fermion number (to — n), and 
according to the Lepage-Savage analysis, we expect 



cx e 



\m—n\M-\-(m-\-n)m 7r /2]T 



(33) 



We will see that this expectation is born out in the A/V formulation, which is the one with a QCD-like sign problem 
at nonzero chemical potential. 

We should not compute graphs contributing to the valence fermion self energy, but rather use the nonperturbative 
solution to the Schwinger-Dyson equation, replacing the mass term (—to) in the fermion propagator by M from 
eq. (22). It is convenient to have this propagator in a mixed {i,p} representation: 



<p',f 



1 



> + M 



|p,^(27r) 2 <5 2 (p'-p)G(p,<'-;) , 



(34) 



with 



G(p,f) 



2ir —iuj^i — ip ■ 7 + M 



2io r , 



(35) 



where lo p = \J\p\ 2 + M 2 . Note that 



G(0,t) 



-M\t\ ( 1 + efthi 



(36) 



is proportional to a projection operator. Since the dynamical mass M includes all nonperturbative contributions to 
the fermion self energy, it follows that the first cumulant for Xr is just 



Ki.o = Tr 



TG(0, 



— M 



T = (X r 



(37) 
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where z = Tr [r ( "a 71 )] is the wavefunction overlap of our chosen observable with the physical fermion state. Thus a 
Monte Carlo simulation that correctly estimates the value of (Xr) will correctly determine the fermion mass to equal 
M by means of the formula eq. (32), provided that T is chosen so that z is nonvanishing. 

As we are interested in how difficult a Monte Carlo determination of (Ar) might be, we turn next to the variance 
At leading order in 1/N, the sum of diagrams for Ki t i is given by attaching G propagators at zero spatial 
momentum to the legs in the first diagram in Fig. [I] with ends k, I at time t = — r/2 and ends i,j at time r/2, and 
cont rac ting the Dirac indices of each valence fermion line with T. Using our result for the four-point function M. in 
eq. (24) (but dividing by A 2 , since our Ap is the propagator for a canonically normalized fermion), we find that the 
D a part of M. is killed by the Dirac trace and only the part contributes, yielding 



«i,i 



1 

VN 



dt 2 Tr TG(0, r/2 - t 2 ) l5 G(0, t 2 - r/2)rG(0, -r/2 - h) l5 G{0, h + r/2) 



dt 2 e" M (l T - 2t2 l+l T + 24l l)Trr 



l + e(f-t 2 ) 7l 



l + e(*i + |)7i 



dw e ~Mt2-ti) 
2n LKjui) 



G-n{t 2 — ti) 



(38) 



where 



~V~ 



1 

VN 



dw e- iu, (*) 



2tt Ar(w) 2m 7r VN 



-m w \t\ 



(39) 



wit h D n , given in eq. (12) and eq. (14) respectively and we approximated D 7r 1 (w) by the pion pole contribution, 
eq. (14), ignoring the branch cut at k 2 < —AM 2 . If we choose T = 1, we find 



«u,i 



dh r dt 2 e- 2M ^- t ^ + ^Ul + e(--t 2 )e(t 1 + - 
2m„VNJ_ 00 1 \ \2 V V 1 2 



-m„|t 2 — ti| 



2?r 



m w fVN 



(40) 



(41) 



where we have (i) assumed we are near the chiral limit with <C /, M and (ii) assumed we are interested in late 
time behavior, r 3> 1/M. Note that near the chiral limit we are finding Ki,o/ \Jki,i oc exp[(— M + m 7r /2)r] < 1 at 
late time, indicating a severe signal to noise problem, in agreement with the Lepage analysis. 

It is interesting to note that if instead we take V = (1 ± 7 i)/2 then Ki t i vanishes at this order in 1/N. This choice 
of r kills the pion contribution to the variance, and so we would expect a noise-free measurement in this case. We do 
not expect this to persist at subleading order in 1/N, nor do we expect that in real QCD one can decouple baryon 
observables from pions so easily, but it seems worth exploring whether correlating initial and final Dirac indices of 
baryon operators (as with this trace with (1 ± 7 i)/2 on valence quark lines) might be able to improve the signal to 
noise problem in real QCD computations. 



2. Power counting for higher K m ,n cumulants 

Higher cumulants can be computed for the A/V formulation using the equivalent <j/tt diagrams as discovered in 
j jIIIB2| and shown in Fig. |J This is not to say that the same diagrams give the cumulants for the A/V and the 



a /it theories. In Fig. [3] we show the diagrams for several of the lower cumulants, where solid lines are the fermion 
propagators l/{—if + M), and dashed lines are the meson propagators G n and G a (both mesons contributing in 
general) with couplings 1 or respectively, at the vertices. Again, we take all incoming and outgoing 2-momenta to 
be zero, and there is one Yni contracted with each pair of like indices in the graph. As we have chosen a canonically 
normalized fermion propagator with one particular flavor as our observable, there are no factors of N at the meson 
vertices, nor are there any loops giving rise to factors of N . However, each meson propagator costs a factor of 1/N, 
and so n m n cx A/' m+n_1 , since we need a minimum of (m + n — 1) mesons to make a connected graph. Furthermore, 
one can see by cutting the graphs at a fixed time that the minimum mass state that can possibly propagate in a graph 
for K m .n with m > n consists of (m — n) fermions with mass M and n pions. Therefore generically we expect these 
cumulant to scale as 

g— (m— n)Mrg- nm^r 

*m,n ^ 1 jVm+ »_ 1 (TO > Tl) . (42) 

This scaling could be violated if T can be chosen so that the pion does not couple, as discussed in the calculation of 
«i l. Then is replaced by 2M, the mass of a fermion/anti- fermion pair. Note that at late time, the 1/N expansion 
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FIG. 3: Contributions to several low cumulants for the A/V theory to leading order in 1/N, using the o/tt method to compute 
them as developed in QIIIB g| Solid lines represent fermion propagators, dashed lines are mesons. The endpoints are to be 
contracted with Yui , Tjji, T^k'- The lightest intermediate state that can appear in a graph for K m , n with m > n consists of 
(m — n) fermions and n pions, and the sum of their masses determines the r-dependence of the cumulant. 

breaks down in the sense that k,2,o becomes smaller than k 2j 2, for example, so long as one is near enough to the chiral 
limit that < M. This is despite the fact that K2,2 is parametrically smaller by 1/N 2 . However, this breakdown of 
the 1/N expansion will not lead to qualitatively different results because contributions to K m ,n which are subleading 
in TV-counting will not lead to lighter intermediate states than the leading calculation, unless there is some fortuitous 
exclusion of the pion at leading order due to the choice of T that does not persisit at higher order. 

The above scaling implies that the distribution for the real part of the fermion propagator near the chiral limit 
becomes highly symmetric about zero at late time because odd moments (for which (m — n) =/= 0) are seen to fall 
off much more quickly than even moments. This is completely consistent with the Lepage-Savage picture for baryon 
propagator distributions in QCD. As in QCD, it will be very difficult to use Monte Carlo methods for the A/V 
formulation to determine the ground state energy with fixed large fermion number. This is not surprising because, 
like QCD, this theory also has a severe sign problem in the grand canonical formulation for studying states with 
nonzero fermion density. 



B. Noise distribution in the a /it formulation and long tailed distributions 

1. A graphical expansion for cumulants of Yr(a,n) 



We now turn to the task of computing the cumulants for the log of the fermion correlator, Y"r in eq. (31a I in the 
a/n formulation. In this case the cumulants are given by the connected graphs derived from the action 

S Y = NS[a, tt] - isY T [a, tt] (43) 

where S[a, n] is given in eq. The Accounting in this formulation is quite different than for the A/V case since 
the o~ and tt mesons are singlets under the U(N) symmetry rather than N x N matrices. It is convenient to associate 
the expansion of NS[a, tt] in meson fields as vertices labeled by black dots, Fig. [4] There are no black tadpoles, since 
we have solved for the chiral symmetry breaking vacuum; furthermore there are no explicit black two-point vertices 
as these are accounted for by using the full meson propagators. The expansion of isYr[a, tt] in powers of the meson 
fields is represented as white dots, the k th term in the expansion drawn as a white vertex connecting k meson lines; 
white vertices occur with any number of meson lines, starting with zero, and each is associated with a power of is. 
The n th cumulant K n is then given by n! times the sum of connected graphs with n insertions of white vertices, 



since each white vertex brings a factor of (is) and we have the expansion eq. (27). Expanding these graphs in powers 
of 1/N is simple: each black vertex entails a power of N, while each meson propagator gives a factor of 1/N. Loops 
do not give factors of N since the a /tt mesons do not carry U(N) flavor quantum numbers. White dots also do not 
give factors of N. Thus a contribution to K n arising from a graph with n white vertices, b black vertices, p meson 
propagators and £ loops is of order N b ~ p . Since every such graph satisfies p — (b + n) = £—1, we can rewrite the order 
of the graph as N~( n ~ 1+f> . It follows that the leading contribution to K n must come from the sum of tree diagrams 
(£ = 0) with n white vertices. The branches of these connected tree diagrams must end on white tadpoles, and are 
quite limited in number; the first few leading diagrams are shown in Fig. [5j 

The sum of tree diagrams can be regarded as the solution to a classical theory, which in the present case is nothing 
other than the statement that the 1/N expansion is equivalent to solving for the generator of cumulants, In $y(s), as 
a saddlepoint approximation dominated by the classical solution that minimizes the action Sy (see Appendix for more 
details). Before tackling this calculation it is worthwhile to note several simplifying features of these tree diagrams: 
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A- 





FIG. 4: Black vertices arise in the expansion of the fermion loop, NS(a, ir) of eq. |5|) in external meson fields. The tadpole 
was eliminated by vacuum minimization, and the two-point function gives the exact meson propagators, so these vertices start 
with the three-point function. 

i. By choosing a T that is neither 75 nor 7^,75, only the a meson can couple to the white tadpole. Since parity 
implies that vertices conserve pion number mod 2, and all tree diagrams end on white tadpole vertices (Fig. |5j, 
it follows that the only mesons propagating in these tree graphs are a mesons. 



ii. Since we are defining our observable in eq. (31a) to be the log of the propagator of a quark with initial and final 
2-momentum p = 0, the meson at the white tadpole vertex must have p = flowing through it as well. Then, 
because all vertices conserve momentum and the leading graphs in Fig. [5] are all tree diagrams whose branches 
end on white tadpoles, it follows that all of the internal meson lines in the graph must be at p = 0, with only 
nonzero energy flowing through the lines. 

iii. As we show below, the large r behavior of the cumulants arises from graphs with zero energy flowing through 
the white tadpole; thus, due to conservation of three momentum at all vertices, the asymptotic r behavior of the 
cumulants is given by the graphs with 3-momentum p vanishing everywhere within the graph. 

We now demonstrate the last point (iii): that the white tadpole enforces zero energy flow through the diagram at 
late time. To do this we compute the tadpole, assuming energy k\ and two momentum k flowing out of the meson 
line: 



tadpole 



i JZz dr' Tr TG(0, r' + r/2)fo(0, r')G(0, r/2 - r') 



1 

V 

V J Cln) 



^(27r) 2 ,5 2 (0)TrrG(O,r/2) 
12 dT ' t |^(2^ 2 (k) e -^'Mfc) 



-r/2 



^(2.) 2 ^(k) 2sin(fclT/2) ^(fc) 
ki 



d 3 k 
(2tt) 3 



(2ir) 3 5 3 {k)6a(k) 



where G(0,t) is given in eq. (36), provided that 



z = Tr 



l + 7i 



(44) 



(45) 



The S(ki) factor justifies computing the tree graphs in Fig. [5] with zero 3-momentum flowing through it, as long as 
we are only interested in the r — > 00 behavior of the distribution of propagators. However, note that at finite r the 
factor sin(/ciT/2)/fci acts as a filter which still allows k\ < 1/t, or equivalently, which still allows the a field to be 
time dependent on time scales > r. This will be relevant to the next section. 



K X : O 



K 2 : OO 



^3 : OOO + 



^4 : OOOO + 





FIG. 5: Leading contributions in the 1/N expansion of k„ for n = 1, .. .,4 in the a/n formulation. Black vertices are given in 
Fig. 01 white vertices are determined by the expansion of isY [a, it] . Lines represent the exact meson propagators G CTl7r derived 
from eq. 
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2. The generator of cumulants from mean field theory 



As mentioned previously, the sum of tree graphs is nothing other than the effective action Sy — NS—isYr evaluated 
at its minimum, a classical solution for the meson fields (see appendix for further discussion). We have shown that 
for an appropriate choice of T, this solution will in general have vanishing pion field and a spatially constant, but 
time dependent, a(t) field. Finally, we have also shown that the large-r behavior is given by an even simpler solution, 
with a a field that is constant in the two spatial dimensions and over time < r. 

It is straightforward to compute Sy(a) for a constant a field. Note that a is defined relative to its vacuum value 
(er) = (M + m), where M is the actual "constituent" mass of the fermion, while m is the "current" mass appearing 
in the Lagrangian. For the NS{a) part of Sy we need only take the expression eq. ^ evaluated at [a + (o - )), up to 
an overall additive constant: 

NS(a) = NVT ( J- ±j- 3 J —L j (46) 

where V is the spatial volume and T is the time extent of the box. The second part of Sy is given by 



isYp = is In 



' /W e— TrF- ' 



-is((M + a)T-kiz) . (47) 



2ir —iujji + M + a 

The equation for the minimum of Sy is therefore given by 

[NVT a{a + 2M - f) + wrsr]^ = (48) 

with solution 

a Q = X - ((/ - 2M) + y/(f - 2Mp - 4i:ist/(NVT^ (49) 

where we choose the solution that vanishes as s — > 0, to recover the correct chiral symmetry breaking vacuum. 
Plugging this solution back into the effective action, we obtain the generator of cumulants of Yy , to leading order in 
the 1/N expansion and at late time: 

ln$( S ) = -Sy(a ) = is f i + (2M - f)r ~ 1 + (1 ~ ^ sC)V2 



with 



12C 

(50) 



Inz-Mr. (51) 



* NVT(2M - /) 2 ' 
Expanding ln$(s) in powers of s yields the cumulants n n for Yy: 

ki=H, * n>2 = (2 . ( " ~^. )! (2M - f)r C- 1 (52) 
_ (n — 2)\ 

The volume factor in the denominator of £ is easy to understand, arising from the normalization of our one particle 
states; the factor of T is puzzling though, arising from our assumptions of a mean field solution which is constant over 
the entire spacetime volume. This docs not make sense, since there is no need for a to adjust from its chiral symmetry 



breaking minimum long before or long after the correlator has acted. As pointed out in the discussion below eq. (451, 
we should not expect a to be constant over time scales > r, and in fact we should expect the mean field to relax to 
its vacuum value for t < — r/2 and t > r/2. Therefore the factor of T in £ — the temporal size of the box — should 
be replaced by ~ r, the time scale of the correlator, and we have 

7T 

C ~ NV(2M - /)2 (53) 

which is independent of r. Perhaps one can compute a more accurate, time-dependent mean field solution, but we do 
not pursue that here. In an analogous calculation of the Polyakov loop distribution at finite temperature, a strictly 
time independent mean field solution is probably exact in the large N limit, due to the homogeneity in time of the 
operator being measured. 
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We conclude this section by remarking that in the particular limit 

7" 

N ^ oo , r ->• oo , — = fixed (54) 

all cumulants n n vanish for n > 3 and Yy assumes a normal distribution, which is to say, a Monte Carlo simulation 
of the fermion propagator will be sampling a log-normal distribution. With a 2 growing linearly with time, this 
distribution will be very heavy-tailed, with the mode being far from the mean. Standard simulation methods would 
fail for such a distribution, but one could use a cumulant expansion of the Monte Carlo data to obtain an accurate 
measure of the fermion mass, as described in Ref. [3]. 



V. DISCUSSION 



Our motivation for returning to the well-worn Nambu Jona-Lasinio model was to elucidate connections between 
chiral symmetry breaking and the sign problem in lattice QCD at finite chemical potential, without having to deal 
with the complications of asymptotic freedom and confinement. What is particularly attractive about the large-iV, 
three dimensional version of the theory we have studied is that (i) it is analytically tractable to compute features of 
the probability distribution of a fermion correlator, and (ii) the theory has two equivalent formulations, one without 
a sign problem, and one with a QCD-like, exponentially severe sign problem. 

We found that in the QCD-like U A/V" formulation, the fermion determinant was complex and that a Splittorff- 
Verbaarschot argument [H [9] could be made to show that the phase of the fermion determinant had to fluctuate 
wildly for /i > with an expectation value exponentially small in the spatial volume. When looking at fermion 

correlators, the distribution evolves to have an exponentially small mean relative to its width, implying a severe 
signal/noise ratio when sampling the correlator using Monte Carlo methods. Furthermore, the severity of the problem 
is controlled by the difference between the fermion constituent mass M (playing the role of the baryon mass in QCD) 
and the much lighter pion mass m n . This followed the Lepage-Savage scaling argument which has even cumulants 
of the distribution diminishing as a power of exp(— m^r), while the odd cumulants - including the mean - fall off 
proportional to exp(— Mr). It is interesting that in three dimensions one can choose an observable for measuring the 
fermion mass (by a particular choice of the matrix T in eq. (31b)) which eliminates the coupling of the fermion/anti- 
fermion pair to the pion, and thereby eliminates the problem of noise in measurement of the fermion mass. Such a 
trick might be a useful way to reduce the noise in simulations of QCD in four dimensions, even if it cannot eliminate 
it. 

In contrast the "cr/V" formulation with even N has no sign problem at nonzero /i, and the correlator distribution 
— the cumulants of which could be analytically computed — was heavy-tailed, and in a certain limit, log- normal. 
This overlap problem would pose challenges to Monte Carlo sampling if the tail became too long, but this sort of 
problem seems to be less severe than the exponential fall-off of the signal/noise ratio in the A/V formulation as seen 
with the cumulant expansion analysis of Ref. [21 HHJ HH] ■ R should be noted that such distributions have been seen 
in QCD for intermediate times, before any asymptotic pion noise sets in, and has been hypothesized to be related to 
clastic scattering between particles |21j : the volume factors in our expressions for the cumulants, eqs. give 
support to this picture. 

Our analysis should make it obvious that the sign problem encountered in QCD at nonzero chemical potential is 
not a fermion problem, but instead a consequence of interactions. In particular, if the particles being studied can 
exist in a tightly bound state of valence fermions, there is going to be a sign problem — a generic feature of a theory 
with dynamical chiral symmetry breaking, in which a light composite pion emerges as a Goldstone boson. This is 
what happens in the A/V formulation of the NJL model studied here: the fact that the A and V fields will bind a 
fermion/anti- fermion pair into a light or massless pion implies that studies of the fermion correlator will be noisy, and 
that at nonzero chemical potential for the fermion there will be a sign problem. In the a/ir formulation without a 
sign problem, the pion exists as a fundamental field and not as a bound state. 

The lessons learned from this model raise the question: is it possible to introduce a mean field for the pion into 
our formulation of lattice QCD (without changing the theory) so that the pion does not appear as a boundstate of a 
valence quark/anti-quark pair, QQ1 For example, one might add and subtract a four-fermion interaction to the QCD 
Lagrangian; the attractive one could be introduced by a and 7r auxiliary fields, while the repulsive interaction could 
be derived by means of auxiliary A and V fields. Then a valence QQ pair would feel the usual gluon attraction, but 
A/V repulsion, so that they do not bind to form a light pion. Instead, pions would appear as fundamental fields which 
could be created through QQ annihilation, but would not appear as a boundstate of valence quarks. It is expected 
that the conventional sign problem could be ameliorated in such a theory - but this example probably introduces 
other sign problems, a cure perhaps as devastating as the disease. 
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Nevertheless, we believe that inventing a way to introduce the pion into QCD as a fundamental field could be an 
important step toward solving the QCD sign problem and beginning to study the properties of ordinary and dense 
matter from first principles. 



Appendix A: Connection Between Mean Field Calculation and Tree Level Cumulant Diagrams 



In Sec. |IV B 2[ we asserted that the sum of tree graphs contributing to the cumulants of Yr is equal to the effective 
action Sy, evaluated at its minimum. Here, we will show this equivalence in more detail. 

We have two different representations for the generating function, W(s). The first is in terms of a functional 
integral: 



Z = e ~W{s) = J [d(j)] e -NS W +isY W 

The second representation is as the generating function for cumulants: 



W(s) = constant — 



(Al) 



(A2) 



Changing variables, s — Nr, we have 



-W(r) 



[d(j>] e 



-NA(r,4>) 



with A(r, <j>) = S{4>) + irY{4>), and 



EfirY 



(A3) 



(A4) 



We now compute W in a large- iV expansion, which is equivalent to a mean field expansion: 



A^ 



(A5) 



where we have defined 



A^ 



n) 



S n A 



A™ = , 



4>=4>o 



with cf>Q being the classical solution that minimizes A and 



This allows us to write 



-W(r) 



-NAM 



constant -ma(°) 
Vdet NAW 



.7=0 



j=o 



(A6) 



(A7) 



In order to determine the leading contributions to K n , we must locate the terms that are both leading in N and 
r-dependent. 



The constant factor in eq. ( A7 ) is independent of r and does not contribute to n ri 



• The factor in brackets in eq. ( A7| is the sum of connected diagrams whose propagators scale as 1/N and vertices 
as N. Thus, these diagrams scale as 



N V-P = N l-L 



(A8) 



using the topological invariant L + V — P = 1, where V, P, L are the numbers of vertices, propagators and loops 
respectively. These diagrams do not contain tadpoles, as there are no terms linear in <f>. Therefore, this quantity 
only contributes to W(r) at one loop and higher and is thus O(N ) and subleading. 
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• The determinant factor in eq. (A7) contributes to W(r) a term ^Tr lniV^ 2 ) = constant + |Tr ln^ 2 ). The 
constant is A^-dependent, but r-mdependent, whereas the second term is r-dependent, but ^-independent. 
Therefore, the determinant factor comes in at O(N ) and is also subleading. 

• We are left with the e - NA(0) term: this gives NA°, the classical action at its minimum, as the leading contri- 
bution to W, at O(N). 

We will now relate A^°\ the action at the classical minimum, to the diagrams in Fig. [H] From above, we see that 
we can expand W(r) in powers of iV: 

W{r) = const + N [W (r) + N^W^r) + N- 2 W 2 (r) + .. .] , (A9) 
where Wo(r) ~ NA^. We can also write W(r) in terms of cumulants, 

W(r) = const - N V ¥-L (N™' 1 ^) . (A10) 



n=l 



From this, we see that the K n may be expanded as 



1 00 i- 

, (All) 

n-1 Z-^i MP y ' 



Kn N n-1 N P ' 

p=0 



with 



00 / . % n 



m 

n=l 



We have previously identified the sum of tree graphs with n white vertices as the leading contribution to the nth 
cumulant, k n $. Thus, we see explicitly that 

z — ' TV. 

n=l 

which proves our assertion that the sum of tree graphs in Fig. [5] is equal to the effective action Sy — NS(4>) — isY{<j>), 
evaluated at its minimum. 
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